Downloading and cropping EJScreen data

Load the required R packages (install first as needed, tbeptools installation instructions are here) .

library(tbeptools)
library(sf)

Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 8.2.1; sf_use_s2() is TRUE

library(mapview)

Download the relevant file from EJScreen. All files available at https://gaftp.epa.gov/EJSCREEN/2022/?C=S;O=A. The file is downloaded to a temporary directory.

# url with zip gdb to download
urlin <- 'https://gaftp.epa.gov/EJSCREEN/2022/EJSCREEN_2022_Supplemental_with_AS_CNMI_GU_VI_Tracts.gdb.zip'

# download file
tmp1 <- tempfile(fileext = "zip")
download.file(url = urlin, destfile = tmp1, method = "libcurl", mode = "wb")

Unzip the geodatabase that was downloaded to a second temporary directory.

# unzip file
tmp2 <- tempfile()
utils::unzip(tmp1, exdir = tmp2, overwrite = TRUE)

Read the polygon layer from the geodatabase.

# get the layers from the gdb
gdbpth <- list.files(tmp2, pattern = '\\.gdb$', full.names = T)
lyrs <- st_layers(gdbpth)$name

# read the layer
dat <- st_read(dsn = gdbpth, lyrs)

Reading layer `EJSCREEN_Full_with_AS_CNMI_GU_VI' from data source 
  `/tmp/RtmpoCsEkM/fileb5ed2b997103/EJSCREEN_2022_Supplemental_with_AS_CNMI_GU_VI_Tracts.gdb' 
  using driver `OpenFileGDB'
Simple feature collection with 86000 features and 158 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -19951910 ymin: -1617130 xmax: 16259860 ymax: 11554350
Projected CRS: WGS 84 / Pseudo-Mercator

Intersect the layer with the Tampa Bay watershed.

# intersect the layer with the tb watershed
dattb <- dat %>% 
  st_transform(crs = st_crs(tbshed)) %>% 
  st_make_valid() %>% 
  st_intersection(tbshed)

old-style crs object detected; please recreate object with a recent sf::st_crs()
old-style crs object detected; please recreate object with a recent sf::st_crs()

Warning: attribute variables are assumed to be spatially constant throughout
all geometries

View the data using mapview.

# view the layer
mapview(dattb)

The layer can be saved as an RData object if needed. Size is very small (~1mb).

# save the layer as an RData object (~1mb)
save(dattb, file = 'data/dattb.RData')

Unlink the temporary files so they’re deleted when you’re done.

# remove temp files
unlink(tmp1, recursive = TRUE)
unlink(tmp2, recursive = TRUE)